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ABSTRACT 


Ihis thesis deals with the enhancement of video images 
degraded by turbid water viewing conditions. An algorithm by Peli 
and Lim has been used with some success for enhancement, but it 
was found to accentuate noise. The thesis examines a combination 
of the Peli and Lim algorithm with three approaches to 
enhancement. 

First, a Short Space Spectral Subtraction algorithm which 
performs the restoration in the density domain, using an estimate 
for the power spectrum of the given data set. The degraded image 
is divided into many subimages and each subimage is restored 
separately and then combined. 

Next, an algorithm for Image Enhancement and Noise Filtering 
by Use of Local Statistics, which uses the assumption that the 
sample mean and variance of a pixel is equal to the local mean 
and variance of all pixels within a fixed range surrounding it. 

Finally, a median filter for noise reduction , where a given 
pixel of a degraded image is replaced by the median of the pixel 
values in a window surrounding it. 

Combination of the algorithms are applied to degraded images, 
and the results are compared and discussed, in each case. It was 


found that noise smoothing can be achieved with the spectral 


subtraction algorithm, and that the local statistics technique 


ylelded very good contrast enhancement. 
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I. INTRODUCTION 


In this thesis, some algorithms for the enhancement of video 
images degraded by  turbid water viewing conditions are 
implemented. The images used to illustrate the processes, were 
recorded during the recovery procedure in a torpedo testing range. 

The recovery operation is carried on undersea with the aid of 
digging equipment and underwater video cameras equipped with 
strong artificial lights. The equipment is controlled and the 
operation is monitored on video monitors at the surface onboard 
the recovery vessel. The recovery equipment in its attempt to dig 
out the torpedo, stirs up sediment wich visually obscures the 
object of interest and impedes the operation. To date no special 
techniques have been used operationally to process or enhance the 
video image before display. 

When it is of interest to “improve” an image, there are two 
broad types of image manipulation processes that cover all 
operations performed to get such improvement: image restoration 
and image enhancement. 

The goal of image enhancement is to process a degraded image 
so that the result is more suitable than the original image for a 
specific application or aids the human analyst in the extraction 
and interpretation of pictorial information [Ref. 1]. Image 


restoration, on the other hand, is intended to bring a degraded 


image back to an ideal degradation-free image as closely as possible 
[Ref. 2]. This makes the image enhancement problem a 
subjective one; the improvement of the image appearance to the 
human viewer is highly dependent on the viewer himself and on 
the application. What is “good” for one person or application is not 
necessarily “good” for another one. Consequently, the suitability of 
the processed image for a specific application makes image 
enhancement techniques very much problem-oriented. 

The degradation process can result from causes such as, an 
imperfect photographic process, imperfect display devices, poor 
contrast due to environmental conditions, different forms of noise 
(channel, quantization, salt and pepper) and others [Ref. 3]. In 
the same manner, techniques for enhancement are dependent of 
the original degradation process. Processes yielding satisfactory 
results for one type of degradation are not necessarily suitable for 
another. 

The interest of this thesis is on examining different techniques 
to enhance images that are degraded by environmental conditions 
that create a great lack of contrast in the recorded images. This 
work is a continuation of previous research begun on the subject, 
which included the implementation of an adaptive filtering [Ref. 4] 
algorithm for contrast enhancement. This algorithm was applied to 
a set of underwater images recorded in turbid water viewing 
conditions. The images are low contrast degraded by noisy 


background. As will be explained in the next chapter, the 


algorithm yielded very good contrast enhancement but, seemed to 
accentuate the noise. The purpose of this work is to investigate 
some noise smoothing techniques to Бе applied to the 
contrast-enhanced images, an alternative contrast manipulation 
scheme, and to compare the results obtained to those found with 
the algorithm already implemented. Besides the images included in 
the experimental results of this thesis, many more images 
resulting from other combinations and variations of the processes 
and different settings of the parameters were obtained. We have 
included here only those thought to be most representative of the 


results obtained. 


H. OVERVIEW OF PREVIOUS WORK 


A. DESCRIPTION OF THE ADAPTIVE FILTERING ALGORITHM 
An adaptive filtering algorithm for image enhancement by 
Peli and Lim [Ref. 4] has been used for contrast enhancement of 
images. The images are typifyed by a large dynamic range, and 
recorded over a medium, characterized by a much smaller 
dynamic range, which causes those regions of the images with 
very high or very low luminance to be poorly represented. This 
algorithm modifies the local luminance mean of an image and 
controls the local contrast as a function of the local luminance 


mean of the image. 


g(ni,n2) 


LOW-PASS 
FILTERING | f, (n1,n2) LINEARITY 


Figure 2.1. Adaptive Filtering for Image Enhancement 





In Figure 2.1 a block diagram of the Adaptive Filtering 


algorithm is shown . In the figure, f(n,,n,) denotes the 
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unprocessed digital image, and fj(n,,n5), which denotes the local 
luminance mean of f(n,,nj), is obtained by low-pass filtering 
f(n, n2). The sequence fy(n;,n2) which denotes the local contrast 
is obtained by subtracting { (пі, пә) from f(n,,n;). The local 
contrast is modified by multiplying fy(n;,n2) with k(f,), a scalar 
which is a function of f,(ny,n2). The modified contrast is denoted 
by f'y(ny,n2). The specific functional form of k(f,) depends on the 
particular application under consideration; a value of k(fji) > 1 
represents local contrast increase while k(f,) « 1 represents local 
contrast decrease. The local luminance mean is modified by a point 
non-linearity and the modified local luminance mean is denoted by 
fi (ni n2). The specific non-linearity chosen depends on the 
particular application under consideration, and in most application 
problems the non-linearity is chosen so that the overall dynamic 
range of the resulting image is approximately the same as the 
dynamic range of the recording medium. The modified local 


contrast and local luminance mean are then combined to obtain 


g(n,,n5), the processed image. 


B. DISCUSSION OF RESULTS 

In this section the application of the algorithm described in the 
previous section to the images shown in Figure 2.2 is presented. 
Figure 2.2 shows two low-contrast noise-degraded images of 512 X 
512 pixels with each pixel represented by 8 bits. Note that the 


small dynamic range of the luminance of both images make 


in 





Figure 2.2 Original Degraded Images. 
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difficult the appreciation of subtle details in them. In this case, 11 
is desired to increase the local contrast in the images to bring back 
the lost details. To achieve this, k(f,) is chosen as shown in Figure 
2.3 (a) and (b) for the images of Figure 2.2(a) and (b), 
respectively. For the case of Figure 2.2(a), the function is chosen 
to increase the local contrast for the low luminance regions, 
maintaining it relatively constant for medium ranges of local 
luminance, and again enhance it, for the brighter regions of the 
degraded image. The results are shown in Figure 2.4(a), where a 
relative contrast increase can be observed. This is especially 
evident in the “fish” shape close to the lower right corner of the 
image. Also note the increase in background noise resulting in the 


textured background observed in the contrast-enhanced image. 


40 100 160 255 40 100 160 255 


LOCAL LUMINANCE MEAN LOCAL LUMINANCE MEAN 
(a) (b) 





Figure 2.3. K(f,) Function for Processing Images in Figure 2.2 
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Figure 2.4. Images in Figure 2.2 Processed with Adaptive Filtering 
for Image Enhancement Algorithm. 
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In the case of the image of Figure 2.2(b), it is desired to 
retain the original local contrast in the darker regions, but 
enhance it in the regions corresponding to medium and high local 
luminance means. This explains the particular shape of the k(f,) 
function chosen to process this image. In the results, shown in 
Figure 2.4(b),some details difficult to be distinguished in the 
un-processed image, are now more evident. Among them are the 
particular shape of the object near the center of the image and 
the form of the ocean floor below the left-most object. Again, it 
can be observed that there is an increase in the background noise 


throughout the contrast-enhanced image. 
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ПІ. DESCRIPTION OF THE ALGORITHMS 


In this chapter, the theorethical foundations and procedures 
for implementation of the different algorithms used in the 
noise-filtering and contrast enhancement problem are presented. 
We start with a description of noise-smoothing techniques to be 
used as post-processing after the adaptive filtering algorithm. Then 
an alternative method for the contrast enhancement 1s presented. 
Finally, the use of median filtering for  noise-smoothing is 


described 


A. SHORT SPACE SPECTRAL SUBTRACTION 
This algorithm, developed by J. S. Lim [Ref.5] performs the 
restoration in the frequency domain, using a converging solution 
for the power spectrum of a given data set that does not depend 
on the original estimate. The procedure is based on two basic 
assumptions. First, each part of an image generally differs 
sufficiently from other parts so that the image cannot be modeled 
by a stationary random process. Second, the power spectrum о! 
the restored image is estimated by the spectral subtraction of the 
additive noise spectra from the degraded image. 
1. Spectral Subtraction 
Frequency domain techniques for image enhancement use 


iterative procedures for estimating the image power spectral 
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density, continuing until certain fidelity criteria are met. The 
density so obtained, is dependent on both the initial estimate for 
starting the iterative procedure and the stopping criterion. Some 
methods converge more rapidly than others, and sometimes first 
and second order statistics of the degrading function must be 
known. 


Consider a model of image degradation as presented in 


Figure 3.1. In this model f(n,,no) represents a digital image, 
р(п;, п») represents a linear space-invariant point spread function, 
and d(n,,nj) represents an additive noise component. Thus, 
g(n,,nj) represents a blurred version of f(n,,n;) and y(n,,n5) is 


the degraded image. 


g(n1,n2) 


b(n1, n2) а(п1, п2) 





Figure 3.1. Model of Image Degradation 


An image restoration system is shown in Figure 3.2. It 


can be seen that the goal of a restoration system is to process a 


degraded image y(n,,n;) through a noise reduction system to 


estimate g(n,,no) which in turn is processed by a deblurring 
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system to get f(ni,nj;), an estimate of the original blur-free noise- 


free image. 


NOISE REDUCTION| g(ni,n2) | DEBLURRING 


SYSTEM 





Figure 3.2. Image Restoration System 


The purpose of this algorithm is to provide a model for 
the Noise Reduction System, as presented in the Figure. A general 


model for the restoration filter is given by 


Р (00) P 
H(Q,0) - | —— — | (3.1) 
ща. Р (®.®,) + @ Р.(0 0) 


мпеге Р(01,ш») and Pq(€,,95) represent the power spectra of 
g(n,,n;j) and d(n,,n;) respectively, a and f are constants. If B is 
unity, H(@;,@2) corresponds to the parametric Wiener filter, and if 
a also is unity, it reduces to the standard Wiener filter [Ref. 1,2]. 
Other values for a and P define different filtering techniques, such 
as power spectrum filtering and geometrical mean filtering. One 


commonly used approach in the implementation of H(w ,,w 2) is an 


iterative procedure which begins with an initial estimate of 


P,(w4,w2) and then iteratively estimates g(nj,ng) апа Р;(оә1,0у2) 
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until a converging solution or a desirable performance is achieved. 
This iterative procedure may be computationally undesirable. 

The Spatial Subtraction approach analytically obtains a 
converging solution, by estimating Pg(w,,«2) from g(n,n;) as 
(1/k)-|G(w,, w2)|2 and using a value for B of Y2. The complete 


solution for this problem, is given by 


E 2 
| б(®,.®,) | = | ¥(@,,0,) | - ак (0.0) (3.2) 
3 С(0 0) - 3 7(9 0) (3.3) 
Еог | Y(&,. 2 О.К. (0 0) 


and zero otherwise 


Here G(w,,@)) represents the discrete space Fourier Transform of 
g(ny,n>), Y(@,,@) is the discrete space Fourier Transform of 
y(ni,n>), Pg(@,,@2) represents the power spectrum of d(n,n;), 
the additive noise component of the degradation model, and k is а 


scaling factor that normalizes the power and energy spectral 
densities. The phase of g(ni,n;) is estimated by the phase of 


y(n,,n3) and the transform magnitude of g(n,,n;) is estimated by 


a particular form of spectral subtraction. 
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2. short Space Implementation | 

As mentioned before, the basic assumption is that each 
part of an image f(ni,no) generally differs sufficiently from other 
parts so that it cannot be modelled by a stationary random field. 
Thus a short space implementation seems  appropiate. Тһе 
technique consists of dividing the original image into many 
subimages and restoring each subimage separately. The procedure 
involves the application of a short space window function 
wi (ning) to overlapping portions of the degraded image y(n,,n;) 


so that: 


y (n,,n5)*v, 4(n,, nj) - [g(n,,n,)*d(n,,n,)) * vi 4014: 25) (3.4) 


or, equivalently, 


Y вт а OT Ы Ы) (3.5) 


The noise reduction system is then applied to yi j(ni,n2) 


to recover g, (n4, n;). The overall full-size restored image g(n,,n2) 


Is given by 
ж 2L 
si nl 2,2,9, ауа) (3.6) 
i=0 ј=0 
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Тһе window function is desired to be a smooth function to 
avoid some possible discontinuities that may appear at the 


subimage boundaries in the processed image and must satisfy: 


ж ж 


P 2 v, (npn) = 1 (3.7) 


i=0 3-0 


for all y of interest 


K 2K 3K 


¥,(n,) ж,(п,-1) и 


,n,-2L) 





Figure 3.3. 2-D Separable Triangular Window 


The condition in Equation 3.6 guarantees that an image 
can be reconstructed from its subimages. To ensure the smoothness 
of the window, a 2-D separable triangular window similar to that 
shown in Figure 3.3 will be used. Each window overlaps its 


neighboring window by half the window length in each dimension. 


Zu 


B. CONTRAST ENHANCEMENT BY USE OF LOCAL STATISTICS 

This algorithm presents a computational technique for contrast 
enhancement on a two-dimensional image array based on the local 
mean and variance. The algorithm, as developed by J. Lee [Ref. 
6] uses the basic assumption that the sample mean and variance 
of a pixel is equal to the local mean and variance of all pixels 
within a fixed range surrounding it. 

Let x; ; be the brightness of a pixel (i,j) in a two dimensional 
NXN image. The local mean and variance are calculated over a 
(2n +1)X(2m +1) window. The local mean is defined as 

: ien jem 
ru Ж, К (3.8) 


kzi-n 12sj-m 


and the local variance is defined as 


1 1+п +m 

2 
V = (x m.) 
1,j (2n+1)(2m+1) 24. а к1 1, 


The algorithm is designed such that a pixel xj; will maintain 


its local mean, and yet permit its variance to be modified by a 


constant factor times its original variance. The procedure is defined 


by 


(3.10) 
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where k, the gain, is the ratio of new local standard deviation to 
the original standard deviation. This approach has the 
computational advantage that the local variance vj ; ls not 
required and only the local mean mj, needs to be computed. If 
we have k » 1, the image will be sharpened as if acted upon by a 


high-pass filter. If O S k « 1, the image will be smoothed as if 


passed through a low-pass filter. In the extreme case,one has k 


O and xi j1s equal to its local mean mij j . 





Figure 3.4. Contrast Enhancement by Use of Local Statistics 


Figure 3.4 shows a block diagram of the algorithm. The local 
mean is evaluated over a (2n+1)(2m+1) window and subtracted 


from the original brightness of the pixel, xij. The difference is 
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then multiplied by the standard deviation ratio to enhance the 


local contrast. Finally, this result is added to the local mean 


yielding the reconstructed pixel value xj. 


C. NOISE FILTERING BY USE OF LOCAL STATISTICS 

This algorithm, also developed by J.Lee [Ref. 6] is actually a 
variation of the previous algorithm adapted for the noise filtering 
problem. In this algorithm, the a priori mean (or variance) of 
the estimated image is calculated as the difference between the 
mean (or variance) of the noise corrupted image and the mean 
(or variance) of the noise itself. Let z; j be the degraded pixel xj ;. 
Then the degraded pixel can be modelled as the sum of the 


noise-free pixel and a white noise sequence. That is, 


Z =X +W (3.11) 


where wi j is the white random sequence with EL wi, ; ]-0. Define 


the a priori mean and variance of x; ; to be 
x, 2 E(x, .J= Elz, ‚]= 2 (3.12) 
i.] 1] 122 ij | 
and 


А 


Е 2 = 2 
о у ОЕЕО абы (3 18) 


1,) 


respectively. The estimated pixel value xj ; is computed by 


24 


X = X, + К, (Zz, =x (914) 


E - (3.15) 





Figure 3.5. Noise Filtering by Use of Local Statistics 


Figure 3.5 shows a block diagram for the implementation of 


the algorithm. Starting with the degraded pixel z; ;,the local mean 


is evaluated and subtracted to obtain the local contrast. Then, 


both Ср) ап k,; can be evaluated, and the value of k thus 


obtained is multiplied by the local contrast. This in turn is added 


to the local mean to get the restored pixel x; j. 


One important parameter is the size of the window over 


which the local mean is estimated. If the window is too small, the 


2S 


noise filtering algorithm is not effective. If the window is too large, 


subtle details of the image will be lost in the filtering process. 


D. MEDIAN FILTERING 

Median Filtering is a non-linear signal processing technique 
originally developed by J.Tukey [Ref. 7], that is useful for noise 
suppression in images. In one-dimensional form the median filter 
consists of a sliding window encompassing an odd number of pixels. 
The center pixel in the window is replaced by the median of the 
pixels in the window. The median of a discrete sequence a4,,a5,..., 
a, for M odd is that member of the sequence for which (M-1)/2 
elements are smaller or equal in value, and (M-1)/2 elements are 
larger or equal in value [Ref. 2]. For example, if the ordered 
values of the pixels within a window are 80, 90, 200, 110, 120, 
the center pixel would be replaced by the value 110, which is the 
middle point of the sorted sequence 80,90,110,120, 200. 

The concept of median filter is easily extended to two 
dimensions by utilizing a two-dimensional window of some desired 
shape such as a rectangle or a discrete approximation to a circle. 
A two-dimensional LXL median filter will provide a greater degree 
of noise suppression than sequential horizontal and vertical 
processing with LX1 median filters, but two-dimensional processing 
also results in greater signal suppression. 

The median filter has been found to be more effective than a 


linear filter for smoothing images with spiky noise degradations 
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because of extrema rejection by the median [Ref. 8].Furthermore, 
the median filter preserves monotonic step edges, that is, it does 
not blur sharp edges as a linear low-pass filter would. 

Another interesting aspect of median filtering that has been 
studied [Ref. 9], is that of the convergence of a filtered image to 
what has been called a root signa! . A root is a signal which is 
invariant under filtering by a particular median filter. The 
technique consists in making succesive passes of a noise-degraded 
image thru a median filter until an image corresponding to a root 
signal is achieved. Definition of a root signal follows, where a 
window width of 2N+1 is used.This definition uses the following 


ideas: 


1) A constant neighborhood is a region of at least N+1 consecutive 
identically valued samples. 


2) An edge is an increasing or decreasing sequence of samples 
which is immediately preceded and followed by constant 
neighborhoods. An edge cannot contain any constant 
neighborhood. 


3) An uzmpulse is a sequence of at most N consecutive samples 
whose values are different from those of the two surrounding 
regions; the two surrounding regions are identically valued 
constant neighborhoods. 


A signal is a root if and only if it contains only edges 
alternating with constant neighborhoods [Ref. 10]. The median 
filtering preserves edges and constant neighborhoods but eliminates 
impulses. The technique is based on the idea that passing a noise 


corrupted root signal through the median filter a sufficient number 
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of times will produce another root signal, which is usually “close” 
to the original root signal. In reference [Ref. 10], the concept of 
"closeness" and a measure of the number of runs to be used are 
calculated as a function of the image and window sizes. 

In this thesis, à median filter included in the Spider library 
[Ref. 11] is used, which is the implementation of the algorithm 
suggested by Huang, Yang and Tang [Ref. 12]. It consists 
primarily of an efficient way of updating the histogram each time 
that a given pixel is replaced by the median within a window. The 
algorithm is described briefly for a 3X3 window: 

Step 1: Obtain the histogram in the first window [Figure 3.6 
(a)] and find the median MDN. Next,count the pixels 
with values smaller than the median of the window. 
LTMDN is that number. 

Step 2: Shift the window right by one pixel [Figure 3.6(b)] 
and update the histogram and LTMDN. First decre- 
ment the counted values of the histogram equivalent 
to gray level values (g/(aJ, g/(dJ, g/(g) ) of pixels a, d 
and g. That is: 

нІ8Т| 2/(а/ | - umsT[g/(a) ]-1 
ні8Т| (а) | - usTlg/(d) |-1 
uisTl g/(g) | » usrlg/(gJ ]-1 
LTMDN is updated as follows: 
IF g/(a) < MDN , LTMDN = LTMDN-1 
IF g/(d) < MDN , LTMDN = LTMDN-1 
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IF giz) « MDN , LTMDN = LTMDN-1 
In the same way for pixels c',f', and 1', 
JusT|g/(c) ] » usr[g/(c) ]*1 
ussT[g/£)] - usT[g/(7) ]41 
HIST[g/(7) |] = нзт| 877) |+1 
IF g/(c) « MDN , LTMDN = LTMDN+1 
IF g//£/) « MDN , LTMDN = LTMDN+1 
IF g//) « МОМ, LTMDN = LTMDN+1 





Figure 3.6. Window for Median Filtering 


Step 3: Update MDN in the previous window, and obtain a 
median in the new window. Let ITH=(number of 
pixels in the window) /2 
Ze MD Tha then 
LTMDN = LTMDN - HISTÎMDN] 
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otep 4: 


MDN = MDN-1 
Repeat this until LTMDN < [ТН 1s obtained. 
= f LIMON STANE 

MDN = MDN+1 

LTMDN = LTMDN + HISTÍMDN] 
Repeat this until LTMDN + HIST[MDN+1] > ITH is 
obtained. This MDN is the median for the current 
window and substitute it in the output image. 


End when one line is finished. Go to Step 2. 
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IV. EXPERIMENTAL RESULTS 


In this chapter, experimental results for different combinations 
of the algorithms described are presented and discussed. Basically, 
the adaptive filtering algorithm [Ref. 4] was used for contrast 
enhancement, аз described in Chapter ll. As previously stated, 
this algorithm yielded very good contrast enhancement, but also 
tended to accentuate the noise. Thus, methods for noise reduction, 
such as short space spatial subtraction [Ref. 5], noise filtering by 
use of local statistics [Ref. 6], and median filtering are used, as 
post-processors to the contrast enhancement operation. The results 
are shown in next sections. Another technique for contrast 
manipulation, contrast enhancement by use of local statistics [Ref. 
6| was also used, as an alternative to the adaptive filtering 
algorithm. These results are also given and compared to those 


previously obtained. 


A. PROCESS TYPE 1 

This section shows a combination of the adaptive filtering and 
short space spatial subtraction algorithms. In Figure 2.2 two 
degraded images are shown. The images are characterized by poor 
contrast due to the turbid water viewing conditions. These images 


were processed for contrast enhancement using the adaptive 
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filtering technique, but the noisy background was heavily 
accentuated as shown in Figure 2.4. Then the spectral subtraction 
algorithm was used to smooth the noise resulting in the images 
shown in Figure 4.1. Inspection of these images shows that this 
particular algorithm was succesful in smoothing the background 
noise present in the contrast-enhanced images, and, at the same 
time, most of the information is preserved. Particularly note, in 
Figure 4.1(a), how details in the beam-shaped object in the center 
of the image are brought back, as well as those of the “fish” in 
the lower right corner of the image. The same observation can be 
made about Figure 4.1(b), in which the noise reduction was 
achieved, at the cost of some signal degradation, in this case. 

This algorithm requires the power spectral density of the noise 
to be estimated. In this case, this was accomplished by taking a 
portion of the noisy background, without any objects and using it 
to estimate the spectrum. Figure 4.2 shows a representation of the 
estimated spectrum of the noise density thus obtained. In the 
figure, the magnitude of the power spectrum is represented as the 
intensity of a 2-D image, each of the image dimensions being the 
corresponding spatial frequency coordinates. 

For the short space implementation, it is required to process 
the degraded image by dividing it into subimages and processing 
each subimage separately. In dividing the degraded image, the size 
of the subimage must be such that the image luminance in it can 


be approximated by a stationary random process. In addition, the 
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Figure 4.1. Images in Figure 2.4 Processed with the Short 
Space Spatial Subtraction Algorithm 
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Figure 4.2. Power Spectrum of the Noise in Figure 2.4(a) 


window function must be a smooth function in order to avoid 
possible discontinuities that may appear at the subimage boundaries 
in the processed image. Thus, two types of smooth functions were 
used- separable 2-D triangular and Hanning windows-of size 32X32 
pixels, overlapped with its neighboring window by half the window 
duration in each dimension. It was noted that there was no 
appreciable difference in using either window  -triangular ог 
Hanning. 

The parameter k, appearing in Equation 3.1 is a scaling factor 


that normalizes the power spectral densities, and was evaluated as 


ша 3 > е) (4.1) 


1,=— 1-ш-< 


where w(li,l2) is the window function used. 
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B. PROCESS TYPE 2 


In this section, a combination of the adaptive filtering for 
Image enhancement and noise filtering by use of local statistics 15 
presented. The images in Figure 2.4 were processed with the noise 
fültering algorithm, using different values for the noise variance. 
Based in the data used for generating the noise power spectrum 
shown in Figure 4.1, a noise variance of 60 was estimated, for the 
image in Figure 2.4(a). For the image in Figure 2.4(b), a noise 
variance of 45 was estimated. The algorithm was found to be quite 
sensitive to the size of the window used to evaluate the local 
statistics, the a posteriori mean and variance. Window sizes of 
8X8, 16X16, and 24X24 pixels were used without noticing any 
appreciable improvement in the processed image. Figure 4.3 shows 
the results of processing the images in Figure 2.4 using a window 
size of 32X32 and the noise variances mentioned above. Note that 
a relative reduction in the background noise was achieved,with the 
processed image resembling the effect of low-pass filtering in the 
background, but without blurring the edges of the objects in the 
images. When larger sizes for the window were used, some 


distortion in the images begun to appear. 


C. PROCESS TYPE 3 
In this section, the result of applying a median filter to a 
contrast-enhanced image, using the adaptive filtering algorithm, is 


presented. The image in Figure 2.4(a) was processed with the 
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Figure 4.3. Images in Figure 2.4 Processed with Noise Filtering 
by Use of Local Statistics Algorithm 
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median filter available in the Spider Library, using a window size 
of 5X5, and making succesive iterations of the processed image, 
until appreciable image degradation start to show up. The image 
obtained after 8 iterations is shown in Figure 4.4. It can be 
observed that some noise smoothing was obtained. Even when the 
resulting image is not noise-free, the remaining noise is of much 
lower spatial frequency than that of the image before the 
median-filtering, and thus seems less objectionable. 

Detailed characteristics of the signal are preserved and some 
details have actually been enhanced. This is true in the case of the 
fish in the lower right corner of the image. A comparison with the 
use of 3X3 window for the filter operation led to the conclusion 
that the smaller window size gives slightly more fidelity in the 
signal preservation. However, since convergence for the 3X3 
window is much slower, this choice is undesirable. A reduction in 
the dynamic range of the median-filtered image can also be 
observed, which presents an effect similar to that of low-pass 
filtering in the background region of the image. However, this 
procedure does not have the undesirable characteristic tendency of 


the low-pass process to blur the edges of the signal. 


oy, 
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Figure 4.4. Median Filtering of Image in Figure 2.4(a) 


D. PROCESS TYPE 4 

In this section the application of the contrast enhancement by 
use of local statistics algorithm is presented. In Figures 4.5 and 
4.6 four low-contrast noise-degraded images are shown. All of 
them are 512 X 512 pixels in size and each pixel is represented by 
8 bits. The images in Figures 4.5 (a),(b) and 4.6(a) are 
characterized by a relatively dark background with a diffused 
effect around the objects in the pictures resulting in low contrast. 
The image in Figure 4.6(b) has an almost uniform brightness 
throughout the picture, which makes it difficult to distinguish 
between the fish and beam shapes from the uniform background. 

Figure 4.7(a) corresponds to the image in Figure 4.5(a) 


processed with the contrast enhancement algorithm, using a value 
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of 4 for the gain k of Equation 3.10. (Recall that k is the ratio of 
the local standard deviation of the processed image to the original 
standard deviation). Note the enhancement in contrast achieved, 
which is particularly evident in the cable connecting the two 
objects in the lower part of the image. 

The image in Figure 4.5(b) was processed using a value of 
k=6, with results shown in Figure 4.7(b). The higher value was 
used to try to enhance the fish shape in the lower right part of 
the image. This enhancement was achieved to some extent. 

Figure 4.8(a) shows the result of processing the image in 
Figure 4.6(a). In this case, a value for k of 4 was again used, 
producing an image that is a much crisper version of the 
corresponding degraded image. Also note that details, such as the 
particular shape of the object in the picture, are enhanced. The 
noise in the background is not greater than that in the original, 
which means that the desired result of contrast enhancement was 
achieved without a significant noise increase. 

The results of processing the image in Figure 4.6(b) are shown 
in Figure 4.8(b). In this case a value of k=5 was used, and the 
contrast of the two objects in the picture are increased with 


practically no increase in the noise. 
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Figure 4.5. Original Degraded Images 
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(b) 
Figure 4.6. Original Degraded Images 
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Figure 4.7. 
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(b) 


Images in Figure 4.5 Processed with the Contrast 


Enhancement by Use of Local Statistics Algorithm 
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Figure 4.8. Images in Figure 4.6 Processed with The Contrast 
Enhancement by Use of Local Statistics Algorithm 
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V. CONCLUSIONS 


Restoration of images degraded by turbid water was considered 
in this thesis. Previous work on the subject included the use of an 
adaptive filtering algorithm for enhancing the local contrast. This 
algorithm was applied to images characterized by poor contrast due 
to the turbid water viewing conditions. The technique yielded good 
contrast enhancement, but tended to accentuate the noise. 

In our work some techniques for noise smoothing were used as 
post-processors to the contrast enhancement operation. 

The specific techniques used were short space spectral 
subtraction, noise filtering by use of local statistics, and median 
filtering. The spectral subtraction algorithm yields very good noise 
reduction, but tends to introduce some signal degradation. The 
filtering technique that uses local statistics -local mean and 
variance- produced noise smoothing in the background. This 
technique was found to be very sensitive to the size of the window 
in which the local statistics were evaluated. Median Filtering did 
not eliminate the noise completely, but the remaining nolse was of 
much lower spatial frequency than the original noise and thus 
seemed less objectionable. Better noise filtering was achieved when 
the degraded image was median-filtered iteratively until some 


signal degradation began to appear. 
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As an alternative to the adaptive filtering algorithm, ап 
algorithm for contrast enhancement by use of local statistics was 
used. This particular technique yielded very good results without 
significant noise increase. 

All processes implemented in this thesis proved to be quite 
computationally intensive. This was especially true for the 
techniques that require the evaluation of local statistics. The 
ultimate goal for the problem addressed in this and previous work, 
15 to perform the enhancement in real time. This suggests a 
natural continuation for the work done so far on this topic, in 


order to implement the algorithms for enhancement in real time. 
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APPENDIX 
COMPUTER PROGRAMS 


PROGRAM SPATSUB 


THIS PROGRAM IS THE IMPLEMENTATION OF THE SHORT SPACE 
SPECTRAL SUBTRACTION ALGORITHM 


ВҮТЕ БҮТЕМСІ(5 12,510) 

CHARACTER* 1 COEN 

CHARACTER* 40 INFILENAME,OUTFILENAME 
INTEGER* 4 IMGBUFF(S12,512), INTARRAY(512,512) 


*** INPUT FILENAMES ; INPUT & OUTPUT *** 
TYPE 10 

FORMAT (C INPUT FILENAME => ',$) 

ACCEPT 15, INFILENAME 

FORMAT (A40) 


ACCEPT 67, IWY 

FORMAT (12) 

TYPE 68 

FORMAT (' COLUMN SIZE OF WINDOW => '.$) 
IWX =16 

ACCEPT 69, IWX 

FORMAT (12) 


IERR - 0 

* ** INPUT ORIGINAL IMAGE *** 

CALL INPUTIMG (INFILENAME,BYTEIMG, IMGX, IMGY) 

* * * CHANGE IMAGE FROM BYTE TYPE TO INTEGER TYPE *** 


CALL BYTE. TO INTEGER (BYTEIMG, IMGBUFF, IMGX, IMGY, IERR) 
IF (IERRR .EG. 1) GO TO 100 


*** PERFORM THE SPATIAL SUBSTRACTION *** 
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CALL CORTO (IMGBUFF, INTARRAY, IMGX, IMGY, IWY, IERR) 
IF (IERR .£Q. 1 ) 60 TO 100 

* * * CHANGE IMAGE FROM INTEGER TYPE TO BYTE TYPE *** 
CALL INTEGER_TO_BYTE (INTARRAY, BYTEIMG, IMGX, IMGY,IERR) 
IF (IERR . EQ. 1) GO TO 100 

*** OUTPUT THE PROCESSED IMAGE *** 

CALL OUTPUTING (OUTFILENAME, BYTEIMG, IMGX, IMGY) 
STOP 

CONTINUE 

TYPE *, !!FORCED TO EXIT !!' 

STOP 


END 
SUBROUTINE CORTO (ENTRADA, SALIDA,IMGX,IM6Y.,IWY ERR) 


THIS SUBROUTINE PERFORMS THE SPATIAL SUBTRACTION ALGORITHN 
USING THE SHORT SPACE IMPLEMENTATION TECHNIQUE FOR IMAGE 
RESTORATION 
INTEGER* 4 ENTRADA (IMGX,IMGY), SALIDACIMGX,IMGY), INTE(32,32) 
REAL*4 AR(32,352),Al(32,32),BR(S2,32),BI(32,32), 
ST(64), CT(64), LBR(64), BRR(S2.32), WINDOW(32,32), 
RENTRADA(S 12,912), IMGWIN(32,32), IARRAY(32,32), MAG(32,32), 
PHAS(32,32),PSD(32,32) 
DOI = 1,IMGX 
DO J = 1,IMGY 
RENTRADA (1,J) * FLOATJ(ENTRADA(I,J)) 
END DO 
END DO 
TYPE *, INTEGER IMAGE,CONVERTED TO REAL DATA’ 
CALL MCHECK(RENTRADA ,IMGX) 
*** BUILD UP THE WINDOW FUNCTION * * * 
CALL BLTWINDOW (16,WINDOW) 
TYPE *, WINDOW FUNCTION BUILT 


CALL MCHECK(WINDOW,32) 
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AK = 0.0 


DO | = 1,32 
DO J =1,32 
AK = AK*(WINDOW(I,J)* *2 
END DO 
END DO 


TYPE *,'AK=>... АК 
**® INITIALIZE OUTPUT ARRAY *** 


DO L = 1,IMG 
DO H = 1,IMGY 
SALIDA (L,M) = 0 
END DO 
END DO 


TYPE *,OUTPUT ARRAY INITIALIZED' 

CALL MCHECK(SALIDA,IMGX) 

*** ESTIMATE POWER SPECTRUM DENSITY OF THE NOISE *** 
CALL ESTINOISE(RENTRADA,PSD) 

TYPE *,POWER SPECTRUM DENSITY OF THE NOISE ESTIMATED 
CALL MCHECK(PSD,32) 

*** MAIN PROCESS SECTION *** 

INSERT LIMITS HERE 


DO l= 1,31 
DO J= 1,31 
LL = ((l-1)* 16)}+1 
MM » ((J-1)* 16041 
DO L = 1,32 
DO M = 1,32 
IA =LL+L-1 
B = MM+M-1 
IF(IA .GE. 16 .AND IÀ LE. 112 .AND. 
IB .GE. 16 .AND. IB LE. 112) THEN 
IMGWIN(L,M)*RENTRADA(IA,IB) 
IARRAY(L,M)EWINDOW(L ,M)* IMGWIN(L P1) 


ELSE IF(IA .6Е.16 .АМОЛА LE. 112) THEN 
IARRAY(L,M)e O0.S* RENTRADA(IA,IB) 

ELSE IF(IB .GE. 16 AND.IB LE. 112) THEN 
IARRAY(L,M)*0.5 * RENTRADA(IA,IB) 
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ELSE 
IARRAY (L,M)ÆRENTRADA(IA,IB) 


END IF 
END DO 
END DO 


TYPE *, DATA ARRAY HAS BEEN WINDOWED' 
CALL MCHECK (IARRAY,32) 
*** TRANSFORM DATA ARRAY INTO COMPLEX FORMAT *** 


DO L=1,32 
DO M=1,52 
AI(L,M) = 0. 
AR(L,M) = IARRAY(L,M) 
END DO 
END DO 


TYPE* CHECK REAL PART’ 


CALL MCHECK (AR,32) 


TYPE *,CHECK IMAGINARY PART (ZEROS) 
CALL MCHECK (AI,32) 


*** EVALUATE 2-D FFT OF DATA IMAGE *** 


CALL FFTS2(AR,ALBR,BI,32,52,ST,CT,LBR,64,2,JERR) 
CALL MCHECK (BR,32) 
DO IX = 1,52 
DO IY = 1,32 

BRRUIX,IY=CABS (CMPLX(BRUIX,IY),BICIX,IY))) 
END DO 
END DO 


TYPE *, FFT DONE,CHECK MAGNITUDE’ 

CALL MCHECK (BRR,32) 

CALL DISPECT (BRR) 

** *DERFORFM SPATIAL SUBSTRACTION * * * 
ALFA =1.5 

PSD = 40.0 

DO IX = 1,32 


DO IY = 1,32 
BR(X,IY) 7 (BRROX,IY))**2) - ((ALFA/AK) * ((PSD(IX,IY) **2))) 
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IF (BR(IX,IY) LT. .0.) THEN 
MAG (IX,IY) = 0 
ELSE 
MAG(IX,1Y) = SQRT(BR(IX,1Y)) 
END IF 
END DO 
END DO 


TYPE *, SPATIAL SUBT.DONE CHECK RESULT 
CALL MCHECK (MAG) 
*** EVAL. MAGNITUDE OF INV, FFT *** 


DO IX= 1,32 
ПО ше Бай 
BR(IX,IY) * MAG(X,IY)* COS(PHASUX,IY) 
BI (IX,IY) = MAG(IX,IY)* SIN(PHAS(X,IY) 
END DO 
END DO 


CALL FFTS2(BR,BI,AR,AI,32,32,ST,CT,LBR,64,-2,JERR) 


TYPE* INV. FFT DONE, CHECK REAL PART 
CALL MCHECK(AR) 


TYPE * ‘CHECK IMAGINARY PART (ZEROS) ° 
CALL MCHECK(AI) 


*** BACK TO INTEGER FORMAT ** * 


DO IX = 1,32 
DO IY = 1,32 
INTECX,IY) 7. JINTCARCX,IY )) 
END DO 
END DO 


TYPE*, INTEGER CONVERSION DONE,CHECK IT 
CALL MCHECK (INTE) 


*** BUILD UP THE OUTPUT ARRAY * * * 


DO IX = 1,32 
DO IY = 1.32 
LL =((l-1)* 16} 1 
MM((J-1)* 16} 1 
IXX zLL*IX-1 
IYY =MM+IY-1 
SALIDA (IXX,IYY) = SALIDACIXX,IYY HINTE(IX,IY) 
TYPE*,IXX AND IYY <> ХХ МУ 
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ЕМО 00 
END DO 


END DO !! END DO J 
END DO !! END DO I 


TYPE*,MAIN PROCESS DONE’ 

CALL TRUNCATE(SALIDA,IMGS,IMGY) 

CALL SCALE_ADJUST(SALIDA,IMGX,IMGY,0.0,255.0) 
TYPE *, SPATIAL PROCESS DONE’ 

RETURN 


END 


SUBROUTINE BLTWINDOW(IWY NIN2) 
REAL* 4 NI(32), NIN2(32,32), N2(32,32) 


L=0 
Pl = 4* ATAN (1.) 
DO | * 1,32 
|! = |-1 
Nil) = 055 - (0.5*COS(2*PI*1II/31)) 
М2(1) = NIC) 
ТҮРЕ ”, МІСІ), N2Q) 
IF (I .GT. 16) THEN 
L=L+] 
NI(I) 21 - (2*L) 
ELSE 
Nid) = 1 
END IF 
END DO 
| =0 
DOI= 1,32 
DO J = 1, 32 
IF (J .GT. 16) THEN 
L=L +1 
NIN2CIJ) = NIC) * (J-2* 1) 
ELSE 
NIN2CI,J) = NIG) * J 
END IF 


NIN2(1,J) = N1(1)*N2(J) 
TYPE *, NIN2(1,J) 
END DO 
L=0 


© 


O O O O 
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END DO 
RETURN 


END 


SUBROUTINE ESTINOISE (RENTRADA,PSD) 


THIS ROUTINE ESTIMATES THE POWER SPECTRAL DENSITY OF 

THE NOISE IN AN IMAGE, BY TAKING THE MAGNITUDE (SQUARED) 
OF THE DISCRETE SPACE FOURIER TRANSFORM OF A (NOISY) PART 
OF THE BACKGROUND. 


TO CHOOSE APPROPIATE REGION, SET PARAMETERS | AND J 
ACCORDINGLY. 


REAL*4 RENTRADA(S12,512), PSD(32,32), AR(32,32), Al(32,32), 
BR(32,52), BI(52,52), ST(32,32), CT(S2,32), LBR(64) 


* ** TAKE A PORTION OF THE IMAGE BACKGROUND *** 


4-80 
K=80 
00 1=1,32 
DO М=1,32 
LL=L+J 
MM=M+K 
AR(L,M)RENTRADA(LL ММ) 
END DO 
END DO 


* * * |NITIALIZE IMAGINARY PART OF DATA *** 


DO L-1,52 

00 М-1,52 
Al(L,M)=0. 

END DO 

END DO 


*** EVALUATE MAGNITUDE OF 2-D FFT *** 

CALL. FFTS2(AR,AI,BR,BI,32,52, ST, CT.LBR,64,2 JERR) 
ГО ЕК2 

00 М=1,32 


PSD(L,M)=CABS(CMPLX(BR(L,M),BI(L,M))) 
END DO 
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END DO 
RETURN 


END 


SUBROUTINE INPUTIMG(INFILENAME,BYTEIMG,IMGX,IMGY) 


READ IN INPUT IMAGE DATA IN BYTE TYPE 

ASSUME INPUT DATA FILE IS A DIRECT ACCESS FILE 
BYTE BYTEIMG(IMGY ,IMGX) 

СНАВАСТЕВ * 40 INFILENAME 


OPEN(UNIT=2, FILE=INFILENAME, STATUS=OLD’, ACCESS=DIRECT’ , 
RECORDTYPE=FIXED’) 


DO I21,IMGY 
READ(2'1) (BYTEIMG(I, J), J=1, IMGX) 
END DO 


TYPE *,' * INPUT IMAGE HAS BEEN READ ` 
CLOSE(2) 


RETURN 
END 


SUBROUTINE BYTE. TO. INTEGER (BYTEDATA, INTDATA, IX, IY, IERR) 


THIS PROGRAM CHANGE DATA (USUALLY 2D IMAGE DATA) IN THE BYTE DATA 
TYPE INTO THE INTEGER TYPE 


BYTE ВҮТЕОАТА (ІҮ, ІХ) 
INTEGER* 4 INTDATA (IY, IX) 


001-1,1Ү 
DO J = 1, IX 


IF (BYTEDATA (I,J) .GE. -128 

AND. BYTEDATA (IJ) LT. O) THEN 
INTDATA(I,J) = BYTEDATA (I,J) + 256 
ELSE IF (BYTEDATA(I,J) .GE.O 

AND. BYTEDATA(I,J) LE. 127) THEN 
INTDATA(I,J) - ВҮТЕрАТА(1,/) 

ELSE 


95 
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100 


100 


ТҮРЕ" ,ЯОАТА ООТ OF RANGE! ROW, COL, DATA => ', 
|, J, BYTEDATA(I,J) 

IERR = 1 

50 TO 100 

ENDIF 


END DO 
END DO 


TYPE*, 'BYTE TO INTEGER CONVERSION HAS BEEN DONE 
CONTINUE 


RETURN 
END 


SUBROUTINE INTEGER TO BYTE(INTDATA, BYTEDATA, IX, IY, IERR) 


THIS PROGRAM CHANGE DATA (USUALLY 2D IMAGE DATA) IN THE BYTE DATA 
TYPE INTO THE INTEGER TYPE 


BYTE BYTEDATA(IX, ІҮ) 
INTEGER* 4 INTDATA(IX, IY) 


DO I21,IX 
DO Jz1,IY 


IF (INTDATACI,J) .LT. O .OR.. INTDATA(I,J) .GT. 255) THEN TYPE*,R 
DATA OUTOF RANGE! ROW,COL,DATA => ', 

|, J, INTDATA(I,J) 

IERR = 1 

60 TO 100 

ELSE IF (INTDATA(I,J) .1E. 127) THEN 

BYTEDATA(I,J) = INTDATA(I,J) 

ELSE 

BYTEDATA(I,J) * INTDATA(I,J) - 256 

ENDIF 


END DO 
END DO 


TYPE *, INTEGER TO BYTE CONVERSION HAS BEEN DONE’ 
CONTINUE 


RETURN 
END 


54 


ООО ОО» 


20 


25 


90 


99 


60 


65 


66 


67 


68 


69 


70 


PROGRAM LOC_STAT 


THIS PROGRAM IS THE IMPLEMENTATION OF THE IMAGE ENHANCEMENT AND 
NOISE FILTERING ALGORITHM 


BYTE BYTEIMG (512,512) 

CHARACTER* | COEN 

CHARACTER* 40 INFILENAME, OUTFILENAME 

INTEGER*2 IMGBUFF(512,512), LAVEIMG(S12,512), LOCMEAN(S 12,512), 
INTARRAY(512,512) 


* ** |NPUT FILENAMES : INPUT & OUTPUT *** 
TYPE 10 

FORMAT ( INPUT FILENAME => ,$) 

ACCEPT 15, INFILENAME 

FORMAT (A40) 

TYPE 20 

FORMAT (° OUTPUT FILENAME => ,&) 

ACCEPT 25, OUTFILENAME 

FORMAT (A40) 


*** INPUT ORIGINAL IMAGE SIZE * ** 


TYPE 50 

FORMAT С ROW SIZE OF IMAGE => ',$) 
ACCEPT 55, IMGY 

FORMAT (14) 

TYPE 60 

FORMAT (* COLUMN SIZE OF IMAGE => ,$) 
ACCEPT 65, IMGX 

FORMAT 914) 


* x * INPUT WINDOW SIZE * * * 


TYPE 66 
FORMAT (' ROW SIZE OF WINDOW »» ',$) 
ACCEPT 67, IWY 
FORMAT (12) 
TYPE 68 
FORMAT (' COLUMN SIZE OF WINDOW => ',$) 
ACCEPT 69, IWX 
FORMAT (12) 


* * * SELECT CONTRAST ENHANCEMENT OR. NOISE FILTERING * * * 
TYPE 70 
FORMAT ( ‘CONTRAST/NOISEFIL ? (C/N) => ', $) 


25 
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80 


85 


90 


95 


АССЕРТ 75. СЕ 
FORMAT (A1) 


IERR = 0 

жек {МРОТ ORIGINAL IMAGE *** 

CALL INPUTIMG (INFILENAME, BYTEIMG, IMGX, IMGY) 

* * CHANGE IMAGE FROM BYTE TYPE INTEGER TYPE* * * 


CALL BYTE TO.INTEGER ( BYTEIMG, IMGBUFF, IMGX, IMGY, IERR) 
IF (ERR £Q.1) GOTO 100 


*** CALCULATE LOCAL MEAN ARRAY *** 


CALL LOCAL. ME AN(IMGBUFF, LOCMEAN, IMGX, IMGY, IWX, IWY, IERR) 
IF(IERR .EQ.1 ) GOTO 100 


Е (СЕ Ед. С ОК. CE EQ. c) THEN 


* * * INPUT K = STANDARD DEVIATION RATIO * * * 

TYPE 80 

FORMAT ( STAND. DEV. RATIO K=>’, $) 

ACCEPT 85, K 

FORMAT (14) 

CALL ENHACONT (IMGBUFF, LOCMEAN, INTARRAY, IMGX, IMGY, K, IERR) 
ELSE 


* * МРОТ МУ = NOISE VARIANCE * * * 
TYPE 90 
FORMAT ( NOISE VARIANCE NV => ', $) 


ACCEPT 95 

FORMAT (14) 

CALL NOISE_FILT ( IMGBUFF, LOCMEAN, INTARRAY, NV, IMGS, IMGY, IERR) 
ENDIF 


SCLMIN = 0.0 

SCLMAX = 255.0 

CALL SCALE_ADJUST (INTARRAY, IMGX, IMGY, SCLMAX, SCLMIN) 
CALL TRUNCATE (INTARRAY,IMGX,IMGY ) 

* * * CHANGE IMAGE FROM INTEGER TYPE TO BYTE TYPE *** 


CALL INTEGER TO. BYTE ( INTARRAY, BYTEIMG, IMGX, IMGY, IERR ) 
ІЕ GERR .EG. 1) GOTO 100 


* * *QUTPUT THE PROCESSED IMAGE *** 
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OOD 


CALL OUTPUTIMG ( OUTFILENAME, BYTEIMG, IMGX, IMGY ) 


STOP 

CONTINUE 

TYPE *.' !! FORCED TO EXIT !! 
STOP 


END 


SUBROUTINE ENHACONT (IMGBUFF,LOCMEAN,INTARRAY,IMNGX,IMNGY,K,IERR) 


THIS SUBROUTINE IMPLEMENTS THE CONTRAST ENHANCEMENT BY USE 
OF LOCAL STATISTICS ALGORITHM. 


INTEGER*2 IMGBUFF(IMGY,IMGX) , LOCMEANCIMGY,IMGOX) , 
INTARRAY(IMGY,IMGX) 


*** EVALUATE LOCAL CONTRAST AND MODIFY IT WITH K *** 


DO | = 1 IMGY 

DO J = 1,IMGX 

INTARRAY(I,J) = K*(IMGBUFF(IJ) - LOCMEAN(,J)) 
END DO 

END DO 


*** ADD MODIFIED LOC. CONTRAST AND LOCAL MEAN *** 


DO | = 1,IMGY 

DO J = 1,IMGX 

INTARRAY (I,J) = INTARRAY (1,J) + LOCMEAN (I,J) 
END DO 

END DO 


** * CORRECT ANY OVERFLOW *** 


DO | = 1,IHGY 

DO J = 1,IMGX 

INTARRAY(I,J) 2 K*(IMGBUFF(I,J) - LOCMEAN(I,J) ) 
END DO 

END DO 


* * * ADD MODIFIED LOC. CONTRAST AND LOCAL MEAN *** 
DO | = 1,IHGY 
DO J = 1,IMGX 


INTARRAY (I,J) = INTARRAY (I,J) + LOCHEAN (I,J) 
END DO 
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*** CORRECT ANY OVERFLOW * * * 


SCLMIN = 0.0 
SCLMAX = 255.0 
CALL SCALE_ADJUST (INTARRAY, IMGX, IMGY, SCLMAX, SCLMIN) 


RETURN 
END 


SUBROUTINE NOISE_FILTCIMGBUFF LOCMEAN,INTARRAY,NV,IMGX, 
ІМ6Ү ЛЕВЕ) 


THIS SUBROUTINE IS THE IHPLEHENTATION OF THE NOISE FILTERING 
ALGORITHM, CALLED BY THE NOISE FILTERING AND CONTRAST ENHANCEMENT 
BY USE OF LOCAL STATISTICS PROGRAM. 


REFERENCES: 


INTEGER*2 IMGBUFF(IMGY,IM6X), LOCMEAN(IMGY IMGX), 
INTARRAY(IMGY ,IMGX) 
REAL*8 AK(128,128) 128,128) 


*** EVALUATE THE VARIANCE OF THE IMPUT IMAGE *** 


TYPE *,CHECK IMGBUFF ARRAY 
CALL MCHECK(IMGBUFF ) 


TYPE *,'CHECK LOCHMEAN ARRAY ' 
CALL MCHECK(LOCMEAN) 

DO | = 1,IHGY 

DO J = 1,IMGX 

IDIF = IMGBUFF(I,J) - LOCMEAN(I,J) 
TYPE 16, IDIF 

FORMAT (' DIF 2» ', 110) 

Q(1,J) » CIDIF* *2) - NV 

ТҮРЕ 17, (ЖІ,) 

FORMAT( Q 7» 'F8.2) 

END DO 

END DO 


TYPE *,’ CHECK THE Q ARRAY ' 
CALL MRCHECK(Q) 

* * * EVALUATE THE GAIN AK *** 
DO I= 1,IMGY 

DO J = 1, IMGX 
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DEN » Q(I,J) + NV 
IF ( DEN ЕО. 0.0) ТНЕМ 
DEN = 1.0 
ENDIF 
AK(Q,J) » (Q(I,J) / DEN) 
END DO 
END DO 


*** ESTIMATE THE RESTORED IMAGE * * * 
DO | = 1, ІМбҮ 
DO J = 1, IMGX 
DIF = AK(I,J)*(IMGBUFF(IJ) -(ОСМЕАМ(1,/)) 
INTARRAY(I,J) = INT(LOCMEAN(I,J) + DIF) 
END DO 
END DO 


RETURN 
END 


SUBROUTINE SCALE_ADJUST (IMGBUFF, IMGX, IMGY, SCLMAX. SCLMIN) 
INTEGER* 4 IMGBUFF (IMGX, IMGY) 


MAXVAL = IMGBUFF (1,1) 
MINVAL = IMGBUFF (1,1) 


DO | =1, IMGY 
DO J = 1, IMGX 


IF (IMGBUFF(I,J) .GT. MAXVAL) THEN 
MAXVAL = IMGBUFF(I,J) 


ENDIF 

IF (IMGBUFF(IJ) .1T. MINVAL) THEN 
MINVAL = IMGBUFF(I,J) 

ENDIF 


END DO 
END DO 


TYPE 10, MAXVAL, MINVAL 

FORMAT(MAX &MIN VALUES =>° 218) 

TYPE *,MAXVAL= ', MAXVAL, TMINVAAL* ', MINVAL 
INTVAL = MAXVAL - MINVAL 


SLOPE = (SCLMAX - SCLMIN)/REAL (INTVAL) 
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DO | = 1, IMGY 
DO J = 1, IHGX 


IF (IMGBUFF(I,J) .EQ. HINVAL) THEN 
IMGBUFF (1,J) = SCLMIN 


ELSE IF (IMGBUFF(I,J) .£Q. MAXVAL) THEN 
IMGBUFF (I,J) = SCLMAX 


ELSE 
IDIF = IMGBUFF (I,J) - MINVAL 
IMGBUFF(I,J) » INT(S) OPE*IDIF) * INT(SCLMIN 
ENDIF 
END DO 
END DO 
RETURN 
END 


SUBROUTINE TRUNCATE(IMGBUFF, ІМ6Х, ІМ6Ү) 
TRUNCATE IMAGE GRAY VALUES INTO 0-255 
INTEGER* 4 IMGBUFFCIMGY, IMGX) 


DO | = 1, IMGY 
DO J = 1, IMGX 
IF (IMGBUFF(IJ) LT. 0) IMGBUFF(I,J) = 0 
IF (IMGBUFF(I,J) .GT. 255) IMGBUFF(1,J) » 255 
END DO 
END DO 


RETURN 
END 
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PROGRAM MEDIFILT 


THIS PROGRAM IS THE IMPLEMENTATION OF THE MEDIAN 
FILTERING IN NOISE REDUCTION 


BYTE BYTEIMG(512, 512) 
CHARACTER* 40 INFILENAME, OUTFILENAME 
INTEGER* 4 IMGBUFF(S512,512), IHST(256), INTARRAY(S12, 512) 


** * |NPUT FILENAMES : INPUT & OUTPUT *** 
TYPE 10 

FORMAT (INPUT FILENAME =› `$) 

ACCEPT 15, INFILENAME 

FORMAT (A40) 

TYPE 20 

FORMAT (OUTPUT FILENAME => ',$) 

ACCEPT 25, OUTFILENAME 

FORMAT (A40) 


*** |NPUT ORIGINAL IMAGE SIZE *** 


TYPE 50 

FORMAT (ROW SIZE OF IMAGE => ',$) 
ІМбҮ=5 12 

ACCEPT 55, IMGY 

FORMAT (14) 

TYPE 60 

FORMAT (' COLUMN SIZE OF IMAGE =>’, $) 
IMGX = 512 

ACCEPT 65, IMGX 

FORMAT (14) 


TYPE * ,IMGX,IMGY 

IERR =0 

*** INPUT ORIGINAL IMAGE * * * 

CALL IMPUTIMG (INFILENAME,BYTEIMG, IMGX, IMGY) 

* ** CHANGE IMAGE FROM BYTE TYPE TO INTEGER TYPE *** 
CALL BYTE. TO INTEGER (BYTEIMG, IMGBUFF, IMGX, IMGY, IERR) 


IF (IERR .EQ. 1 ) 60 TO 100 
*** PERFORM THE MEDIAN FILTERING *** 


DO! = 1, 20 
CALL MED! (IMGBUFF, INTARRAY, IMGX, IMGY, 5, 5, IHST, 256) 
DO IX = 1, IMGX 


EM 


100 


DO IY = 1, IMGY 
IMGBUFF(IX,IY) = INTARRAY(IX, IY) 

END DO 

END DO 
END DO 
*** CHANGE IMAGE FROM INTEGER TYPE TO BYTE TYPE *** 
TYPE *, IMGX, IMGY 
CALL INTEGER_TO_BYTE (INTARRAY, BYTEIMG, IMGX, IMGY, IERR) 
ІР ЕКЕ .EG. 1) 60 TO 100 


* ** QUTPUT THE PROCESSED IMAGE *** 

CALL OUTPUTIMG (OUTFILENAME, BYTEIMG, IMGX, IMGY) 
STOP 

CONTINUE 


TYPE *, ! ! FORCED TO EXIT ! !' 
STOP 


END 
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